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. The necessity of accurately taking into account the existence of nonequivalent op- 
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^ ' demonstrated that Bose systems in the presence of the Bose-Einstein condensate and 
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pH , representation for the field operators is constructed, allowing for a self-consistent de- 

O I scription of Bose-condensed systems. Equations of motion are derived from the given 
Hamiltonian, which guarantees the validity of conservation laws and thermodynamic 
self-consistency. At the same time, the particle spectrum, obtained either from diago- 
nalizing this Hamiltonian or from linearizing the field-operator equations of motion, has 

' no gap. The condition of the condensate existence assures the absence of the gap in the 

^ . spectrum, irrespectively to the approximation involved. The suggested self-consistent 
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1 Introduction 



Physics of Bose-condensed systems has received in recent years much consideration, both 
in experiment and theory. Numerous references on current hterature can be found in the 
review articles [1-5] . The history of works on Bose- Einstein condensation is surveyed in Ref . 
[6]. 

The microscopic theory of Bose-condcnsed systems is based on the Bogolubov ideas, 
which have been apphed to a low-temperature dilute gas in the pioneering papers [7,8] and 
expounded in detail in the books [9,10]. A key point in the Bogolubov approach is the 
separation of condensate on the level of the field operators, with replacing the condensate 
operator by a nonoperator quantity, playing the role of a condensate wave function. This 
procedure is asymptotically exact in the thermodynamic limit [11]. In general terms, this 
is formulated as the Bogolubov shift, when the field operator is shifted by the condensate 
wave function. This approach gives for the dilute Bose gas the famous Bogolubov spectrum 
having at low momenta the phonon character [7-10]. Beliaev [12,13] included in the Bo- 
golubov approach interactions between quasiparticles, obtaining, by means of perturbation 
theory at zero temperature, attenuation parameters. Hugenholtz and Pines [14] formulated 
the condition on the chemical potential, providing a gapless quasiparticle spectrum at zero 
temperature. This condition was generalized to finite temperatures by Bogolubov [9,10] and 
Hohenberg and Martin [15]. The Behaev zero-temperature calculations were extended to 
finite temperatures by Popov and Faddeev [16] and Popov [17-19]. The Bogolubov theory 
was applied to inhomogeneous systems by Fetter [20]. The Beliaev-Popov approach was 
extended to nonuniform trapped gases by Fedichev and Shlyapnikov [21] and by Morgan et 
al. [22-25]. 

In the lowest dilute-gas approximations, such as the Bogolubov and Beliaev-Popov ap- 
proximations, the quasiparticle spectrum is gapless, in agreement with the Goldstone and 
Bogolubov theorems [9,10,15]. However, in higher approximations, the spectrum acquires a 
nonphysical gap, which happens, e.g., in the Hartree-Fock-Bogolubov approximation. The 
standard remedy for removing the gap is to construct higher-order vertex equations, such 
as the Bethe-Salpeter equation or Lipmann-Schwinger equation, including additional terms 
into the theory in order to compensate the gap [26]. But so constructed compensations 
are not a result of a self-consistent scheme, in which all equations would follow from one 
generally prescribed procedure. As a consequence, the resulting modifications suffer from 
violating conservations laws and thermodynamic consistency. In order to preserve conser- 
vation laws, all equations in the theory must be derivable from a given Hamiltonian or 
some variational functional [27,28]. For one chosen approximation, it is always possible to 
construct an artificial variational functional that would yield a gapless spectrum and even 
preserving conservation laws [29,30]. However, this cannot be considered as a self-consistent 
theory allowing for a systematic derivation of other approximations. In this way, all the- 
ories of Bose-condensed systems suffer either from the absence of conservation laws and 
thermodynamic consistency or from the presence of an unphysical gap in the quasiparticle 
spectrum. A detailed classification of theories onto either conserving but gapless or gapless 
but nonconserving has been done by Hohenberg and Martin [15]. 

The aim of the present paper is to suggest a self-consistent theory that would be both 
conserving as well as gapless. To be conserving, the theory must be based on one mi- 
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croscopic Hamiltonian, from which all equations would follow. To be gapless, the theory 
requires an accurate mathematical formulation, taking account of the existence of the Bose- 
condensed fraction. Such an accurate mathematical formulation can be done by explicitly 
taking into consideration the appearance of nonequivalent operator representations associ- 
ated with canonical commutation relations [31]. Since the theory to be developed rests on 
the notion of nonequivalent operator representations, it is necessary, first of all. to spec- 
ify its meaning for Bose-condensed systems, which is done in Sections 2 and 3. Then, in 
Section 4, a composite representation is defined, making it possible to accurately describe 
Bose-condensed systems. Equations of motion are derived in Section 5 from a microscopic 
Hamiltonian based on the composite representation. To make the resulting equations more 
transparent, they are illustrated for the case of a uniform dilute gas in Section 6. The gen- 
eralization to nonuniform matter is given in Section 7. Everywhere the system of units is 
employed, where the Planck constant h = 1 and the Boltzmann constant /c^ = 1 are set to 
unity. The material of this paper is based on the Lectures [32] . 



2 Bogolubov shift 

The Bogolubov shift of the field operator is commonly defined as 

V'(r,t) =7/(r,t)+^i(r,t) , (1) 

where ip{r,t) and ■?/'i(r, t) satisfy the Bose commutation relations and i]{r,t) is a condensate 
wave function, assumed to be not an identical zero. In what follows, for the sake of brevity, 
the time dependence often will not be shown exphcitly, so that we shall write -(/^(r), r]{r), 
and ipi{r), though keeping in mind that, in general, these quantities do depend on time. 

Field operators act on functions pertaining to a Fock space. Let the field operators ip{r) 
and ip'^{r) be defined on the corresponding Fock space J-{ip)- Let |0 > be the vacuum state 
of J-{ip) , such that for any r one has 

V'(r)|0>-0. (2) 

Then any function (p e ^('0) can be generated by a repeated action of '0^(r) according to 
the rule 

V=T.^ fn{ri,...,rn)llip\ri) dri\0> . 
n=o vn! j=i 

In that sense, one says that the Fock space J-'{ip) is generated by ip'^i^r). 
It is easy to notice that the state |0 > is not a vacuum for V'i(r), since 

^i(r)|0>= -^(r)|0>^ 0. (3) 

A vacuum for ipi^r) is another state, which will be denoted as |0 >i, and which satisfies the 
vacuum-state equation 

7Ai(r)|0>i=0. (4) 

In line with the general rules [33], the creation operator V'i(r) generates the Fock space 
Clearly, the vacuum state |0 >i is not a vacuum for ip{r), as far as 

^(r)|0>i=r;(r)|0>i7^0. (5) 
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Thus, transformation (1) cannot be understood as an operator equality defined on one Fock 
space. 

Let us introduce the operator 



C' = exp|y r]*{r)ilj{r) — r]{r)^\r) (irj (6) 



and its inverse 



C-^ = exp |- y" [r^*(r)V'(r) - r?(r)V^(r)] drj . (7) 
With the help of these operators, relation (1) takes the form 

^{r) = CM^)C-' (8) 

and, conversely, 

V'i(r) = (7-^(0(7. (9) 

By this definition, the operator C transforms functions from into those in J^{ip), while 

C"^ transforms the states of ^{ijj) into those of In particular, the relation between 

the vacua in J-{ip) and in J-{ipi) is 

|0 >i= (7-^10 > , (10) 

which is obvious from Eq. (9). 

To better understand the connection between the spaces ^(V^) and J-'{ipi), we can use 
the Baker-Hausdorff formula for any two operators A and B, whose commutator [A, B] is 
proportional to a unity operator. Then 

e^+^ = eVexp(- ^ [i, B]^ . 

Using this formula, operator (7) can be represented as 

C~'^ — exp i^J r]{r)ilj^{r) dr^ cxp 1^— J r]*{r)ip{r) dr^ cxp ^ J |?7(r)p(ir| . (11) 

Form here, we have 

(7-^10 > = expj- ^y |77(r)|2cir|exp|y7;(r)V^^(r) (irj |0 > . (12) 

The right-hand side of Eq. (12) is nothing but the coherent state 

|77 > = 770 exp I y r]{r)ip^{r) drj |0 > , 

normalized to unity, so that 

|r;o| = expj- i y \r]{r)\^dr^ . 



4 



By definition, the coherent state is given by the eigenproblem 



^{r)\r] > = r]{r)\r] > , (13) 

whose eigenvalue is a coherent field r]{r). Hence, vacuum (10) is the coherent state, 

|0 >i = (7-^10 > = |r/ > , (14) 

and the condensate wave function in the Bogolubov shift (1) is nothing but the coherent 
field 

77(r) = < r/|V'(r)|?7 > . (15) 
The scalar product of the vacua |0 > and |0 >i is 

< 0|0 >i = < 0|77 > = exp |- ^ J |?7(r)|2 drj , (16) 

as far as 

< 0| exp I J ri{r)ip\r) drj |0 > = 1 . 
The condensate wave function defines the condensate density 

Po(r) = |r/(r)r (17) 
Respectively, the latter determines the number of condensed particles 

No^Jpo{r)dr. (18) 

So, the scalar product (16) gives 

< 0|77 > = exp (- ^ TVo) ■ (19) 

Bose-Einstein condensation implies that the number of condensed particles Nq is large, being 
of the order of the total number of particles N. Therefore product (19) is asymptotically 
zero, 

< 0|r7 > ~ {No > 1) . 

That is, the vacua |0 > and |0 >i are asymptotically orthogonal. Similarly, all states 
generated by V'^(r) from |0 > are orthogonal to the states generated by ipi{r) from |0 >i. This 
is equivalent to saying that the Fock spaces J-'{ip) and J^{ipi) are asymptotically orthogonal 
to each other. There is the sole state that is shared by both these spaces. This is the coherent 
state |?7 > of which is at the same time the vacuum |0 >i= \r) > of J^{ijji). But since 

the spaces J-'{ip) and J-{ipi) are of continuous power, the sole intersection is called to be of 
zero measure. 

The operator transformation (6), from J-'{ipi) to J-'{ip), is not unitary. The inverse trans- 
formation (7), from to is not the same as C'^, since it is not defined in one 
space. To be unitary, an operator and its Hermitian conjugate must be defined in one space. 
Therefore the Bogolubov shift (1) realizes the unitary nonequivalent representations of the 
field operators satisfying the Bose commutation relations. Strictly speaking, Eq. (1) has 
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no sense of an operator equation, since its left-hand and right-hand sides are defined on 
different Fock spaces orthogonal to each other. More correctly, it has to be interpreted not 
as an equation but as a replacement 

^(r)^r?(r)+Vi(r), (20) 

with its left-hand side being defined on J^{ip) and the right-hand side, on ^{ipi). 

As a physical interpretation, we may say that the field operator ipijc) describes the Bose 
system before Bose- Einstein condensation has occurred, while the operator 77(r) -|- ■01 (r) 
corresponds to that Bose system after the occurrence of Bosc-Einstcin condensation. This 
is somewhat analogous to the Van Hove picture of scattering [34,35], according to which the 
particle states, before the scattering process has occurred, are orthogonal to their states after 
the scattering. That is, the Fock space of particles before the scattering is orthogonal to 
their Fock space after the scattering. Thence the particle field operators before and after the 
scattering process pertain to unitary nonequivalent representations of the same commutation 
relations. 

It is important to emphasize that the right-hand side of replacement (20) clearly breaks 
the gauge symmetry. The latter could also be broken in the space J-{il^) by adding to 
the Hamiltonian infinitesimal sources breaking this symmetry. However, irrespectively to 
whether the gauge symmetry in is broken or not, this space remains always orthogonal 

to ^{'ipi). This means that by breaking in J^{ip) the gauge symmetry solely by infinitesimal 
sources (without the Bogolubov shift) one cannot make the transformation from T{ip) to 
For the correct description of a Bose-condensed system it is not sufficient just to 
add to a Hamiltonian infinitesimal sources breaking the gauge symmetry, but it is neces- 
sary to make the Bogolubov shift, transferring the consideration to another space. Recall 
that separating in the Hamiltonian the condensate field operators is equivalent to the Bo- 
golubov shift, because the condensate operators, under the thermodynamic limit, become 
nonoperator quantities representing the condensate wave function [9-11]. 

3 Canonical transformations 

The Bogolubov shift (1) or (20) is an example of a simple canonical transformation re- 
alizing unitary nonequivalent operator representations. Another canonical transformation, 
constantly employed in the theory of Bose systems, is the Bogolubov canonical transforma- 
tion 

Ofc = Ukhk + fcfelfc (21) 

from the operators and a\ to hk and where the label k implies a set of quantum 
numbers, e.g., the momentum. The inverse transformation is 

h = ulttk - vlalf, . (22) 

In order that and bk would satisfy the same Bose commutation relations, the coefficient 
functions have to obey the normalization condition 
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Transformations (21) and (22) also realize unitary nonequivalent operator representations. 
Showing this below, we shall explain how this case differs from that considered in the previous 
section. 

Let |0 >a be the vacuum state for the operator a^, so that 

ak\0 >a= . (23) 

Using the standard procedure [33], one can construct the Fock space T{ak) generated by a\. 
The state |0 >„ is not a vacuum for 6^, since 

&fc|0 >a= -vlalk\Q >a 7^ . 

The operators hk possess another vacuum |0 for which 

fefe|0>b = 0. (24) 
And |0 >b is not a vacuum for a^, because 

ajk|0>6 = 'i;*jfc6Lfe|0>6y^0. 

The Fock space Tihk) is generated by h\. 

Usually, the coefficient functions can be chosen to be real and invariant under the change 
of k onto —A;, which we shall assume below setting 

Introducing the quantity 

ak = ln{uk + Vk) , (25) 

one may write 

Uk — cosha/j , Vk — sinhctfe . (26) 

Define the operator 

-B = exp I ^ J2ak (aka-k - a-fc«I) | (27) 

and its inverse 

B'^ = exp |- ^ ^ Ofc (afcO-fe - alfc4) | . (28) 

Using this, transformations (21) and (22) can be represented as 

ak = BbkB-^ (29) 

and, respectively, as 

bk = B~'akB . (30) 

The operator B transforms T{hk) into J-{ak), while B~^ transforms J-{ak) into J-{bk). Thus, 
the vacua in these spaces are related by the transformation 

\0>h = B-'\Q>a . (31) 
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Involving the Baker-Hausdorff formula, we can represent operator (28) as 
B-^ = exp Q ^afeaLfcolj exp ^ ^afcafeO-fej exp |- ^ Y.^li'^ + 24ait)| , (32) 
keeping in mind its action on vacuum |0 >„, for which 

^akOip a'_kal,apa_p |0 = «fc (l + 24«ik) |0 >a = X! "fcl^ ■ 

kp k k 

Then for vacuum (31), we find 

|0 >6 = S-^|0 >„ = exp (^-\Y1 exp ^ afeO^at j |0 >„ . (33) 
The scalar product of the vacua |0 >„ and |0 is 

„<0|0>,= exp|-- J ^Ij^j , (34) 



where V is the system volume, and the equality 



a < 0|exp 



alj^al ]\0>a= 1 



is taken into account. The integral in Eq. (34) is positive. As is evident, the vacua are 
asymptotically orthogonal, so that 

,<0|0>,,~0 (V^oo). 

The Fock spaces J^{ak) and J^{hk) are mutually orthogonal in the thermodynamic limit. The 
operator B is not unitary. Thus, the canonical transformations (21) and (22) realize unitary 
nonequivalent operator representations. Since the left-hand and right-hand sides of Eq. 
(21) are defined on different Fock spaces, it should be understood not as a straightforward 
operator equality, but as a replacement 

Ofe — > Ukhk + v*_kblk ■ (35) 

The main physical difference between the transformations (20) and (35) is that the latter 
is not related to the existence or absence of Bose-Einstein condensate, while the Bogolubov 
shift (2) transforms a system without Bose condensate to a Bose-condensed system. 



4 Composite representation 

The description of a Bose system without Bose-Einstein condensate can be done in terms of 
one field operator ilj{r) defined on the Fock space generated by ip^{r) from a vacuum 

state |0 >. But a Bose-condensed system requires to employ two variables, the condensate 
wave function r){r) and the operator of noncondensed particles V'i(r) defined on the Fock 
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space ^(V'l) generated by ^l{r) from the vacuum |0 >i. In order to stress that there 
are two variables for a Bose-condensed system, it is possible to treat both condensed and 
noncondensed particles on an equal footing by introducing a composite representation for 
the field operators. 

As is shown above, the condensate corresponds to the coherent state. We may define the 
one-dimensional space 

^o = {\v>} (36) 
containing the coherent state for a field operator ■00 (r), for which 

^o(r)|ry> = vir)\v> ■ (37) 

The operator of noncondensed particles i/Ji{r) is defined on the Fock space J-'{iIji) generated 
by ipi{r) from the vacuum |0 >i. In order to avoid the double counting of the degrees of 
freedom, the orthogonality condition 

J 77*(r)V'i(r)dr = (38) 

is imposed. 

Since in the space J-'{ipi) the gauge symmetry is explicitly broken, the statistical average 
< ipi{r) >T{ti>i) over this space may be not zero. This, however, would result in the noncon- 
servation of quantum numbers. For instance, in the case of particles with nonzero spin, the 
latter would not be conserved. Or, Fourier transforming ipi{r) to Ofc, and having a nonzero 
statistical average < >r(ak)-i one would confront the absence of momentum conservation. 
To eliminate such unpleasant features, one defines the restricted space !Fi C J^{ipi), such 
that the restricted average 

<V'i(r)>^, =0 (39) 
be zero. Formally, the restricted space is defined as 

Ti = {ipe ^(V-i)! < ^i(r) = 0} . (40) 

An explicit and mathematically correct definition of the restricted space (40) can be done 
in terms of the weighted Hilbert spaces [36,37]. 
Let us now introduce the composite Fock space 

T = Tq®Tx, (41) 

being a tensor product of the spaces (36) and (40). And let us define the composite field 
operator 

^P{r)=M^)®M^) (42) 

acting on the composite space (41). In what follows, we shall simplify the notation often 
writing, instead of the sign of the direct summation ©, just the usual signs of summation or 
subtraction. We shall also omit, as is usually done, the unity operators Iq and li defined in 
the corresponding spaces JFq and J-'i. 

The condensate wave function, according to definition (37), is the coherent field 

r]{r) = < ri\i)o{r)\ri > . (43) 
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The number-of-particle operator of condensed particles, defined in the space is 

No = J 4{r)Mr) dr . (44) 
The condensate wave function is normahzed to the number of condensed particles 

No = < r]\No\r] > = j \r]{r)f dr . (45) 

The vectors / e pertaining to the composite space (41), have the structure of the 
tensor product 

f ^\ri>0ip , 

in which \ri >G J-'o and ip E J-'i. For any two vectors fm, fn £ ^, the matrix element of the 
composite field operator (42) is 

/m V^(r)/n ^Vt< V\Hr)\V >Vn^vt H^) + V'lW] ■ (46) 

Respectively, the matrix element of an operator A[iij] from the algebra of local observables 
on !F takes the form 

fm M'4']fn ^^i< Vl^MV >^n = ^i Mv + H¥>n ■ (47) 

The trace of such an operator is 

Tr^i[^] = Tr^, < r/|i[^]|r/ > = Tr^,i[^7 + ^i] . 



The statistical average over the space J^, with respect to a statistical operator p[ip], taken 
at the initial time t = 0, writes as 

< >^ = Ti^mAiTp] = Tr^^p[ri + ^,]A[ri + . (48) 

Defining the average over JF^ as 

< A[r) + iJi] >:f, = Tr^iP[?? + ^M^l + V'l] , (49) 

we get 

<A[i;] >^=<i[r^ + V'i] . (50) 
The operator of the total number of particles 

N = j ^/^t(r)v5(r) dr , (51) 

with the composite field operator (42), keeping in mind the orthogonality condition (38), 
can be represented as the direct sum 

TV = TVo e TVi (52) 

of the condensed-particle operator (44) and the operator for the number of noncondensed 
particles 

Ni = j ^I(r)^i(r) dv . (53) 
10 



The average of the total number-of-particle operator (51), according to Eq. (50), becomes 

< TV >^ = A^o+ < ^1 >^i ■ (54) 

The equahty 

Ar = <iV>^ (55) 

for the total number of particles serves as a normalization condition, additional to the nor- 
malization (45) for the number of condensed particles. Both these normalization conditions, 
(45) and (55), must be taken into account in a self-consistent theory. 



5 Equations of motion 

The energy operator is given by the Hamiltonian 

H = J V5t(r) ^_ ^ + t/j ^(r) dr+^J V5t(r)v5^(r')$(r - r')V^(r')V5(r) drdr' , (56) 

in which, m is particle mass, U — U{r,t) is an external field, and $(r) = $(— r) is an 
interaction potential, assumed to be integrable, so that 



$(r) dr 



< oo 



The field operators ip{r) = ip{r,t) are given in the Heisenberg representation, but the tem- 
poral variable t is not explicitly shown just for making the formulas less cumbersome. 

To take into consideration two normalization conditions, (45) and (55), it is necessary to 
introduce two Lagrange multipliers, /lo and /i, entering the grand Hamiltonian 

H[i^]^H-lioNo-iiN . (57) 

Here n is the conventional chemical potential and /^o is an additional Lagrange multiplier 
guaranteeing the validity of the normalization condition (45). 

With the composite field operator (42), the Hamiltonian (57) can be represented as a 
sum ^ 

n=0 

in which the terms are classified according to the number of the noncondensed-particle oper- 
ators in the products of each term. Thus, the zero-order term contains only the condensate 
field operators, 

//(°) = J 4{v) (-^ + U-^io-^^^ ^o(r) dr+ 

+ i J 4(r)4(r')$(r - r')^o(r')^o(r) drdr' . (59) 

The first-order term is 

//(^) = 1 4(r) ^-^ + U-fx^ V'o(r) dr + l V'i(r) ^-^ + U-fx^ 4(r) dr + 
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+ 1 $(r - r') [V'I(r)V'S(r')Vo(r')V'o(r) + 4{r)^l{r')M^')M^)] drdr' . (60) 

However, a Hamiltonian term of the first order in V'i(r) or t/'i(r), being considered on the 
space J^i defined in Eq. (40), is effectively zero due to the conservation condition (39). If a 
Hamiltonian would contain such linear in ■01 (r) or ■0i(r) terms, condition (39) could not be 
satisfied. The same can be said in another manner. To take into account restriction (39), we 
could add to the grand Hamiltonian (57) one more term with a Lagrange multiplier assuring 
the validity of Eq. (39). To this end, this additional term has to be chosen so that to cancel 
the linear term (60). In any way, 

//(I) = . 



(61) 



For the quadratic term in Eq. (58), we have 

H^'^ = 1 4{r) (-^ + U-iJ^ V'i(r) dr + 

+ 1 $(r-r') [V'i(r)V'I(r')V'i(r')V'o(r) + V'i(r)V'I(r')V'o(r')V'i(r)+ 

+ ^ ijl{r)4{r')Mr')M^) + I 4ir)4{r')M^')Mr)\ drdr' . 
The third-order term is 



<l>fr 



') ^(r)^I(r')^i(r')^i(r)+^I(r)0l(r')0i(r')^o(r)] drdr' . 
In expressions (62) and (63), the symmetry $(r) — •^*(— r) is used. The last term is 

= 1 y ^/;t(r)^/;|(r')$(r - r')^/'i(r>i(r) drdr' . 



(62) 
(63) 

(64) 



In order to preserve all conservation laws, such as the continuity equation, the density 
of current, and the energy-momentum conservation laws, all equations of motion must be 
derived from the same Hamiltonian. A theory with such Hamiltonian-derivable equations 
will automatically be conserving. Since we have two field-operator variables, we should have 
two types of equations of motion. 

The first equation prescribes the evolution of the condensate operator 



I — ipo{r,t) = . 
dt S^ljl{r,t) 



''1 



(65) 
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For the Hamiltonian (58), one has 

64{r,t) 6ijl{r,t) 

In the expressions below, we shall again write tpo{r) and ipiir), omitting the time variable 
in order to simplify the formulas. Then we find the variational derivatives of the zero-order 
term. 



S4{r,t) 



— + [/ - //o - Vo(r) + /$(r - r')4{r')Mr')M^) dr' , 
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of the second-order term, 



$(r - r') Vl(r')V'i(r')V'o(r) + V'I(r')V'o(r')V'i(r) + V5(r')V'i(r')V'i(r) dr' , 



and of the third-order term 



$(r - r')V'I(r')V'i(r')V'i(r) ^r' . 



Substituting this into Eq. (65) and taking its coherent average yields 



(-^ + U-i,o-l^^v{r) + J^r-r')[\ri{r')\'rj{r)+X{r,r')] dr' , (66) 



where 



X{r,r') = 4{r')^^;,{r')rJ{r) + ^^;l{r')rJ{^^^^^ . (67) 

We may note that Eq. (66) can also be derived from the equation 

(68) 



which follows from taking the coherent average of Eq. (65). 

The second equation of motion describes the evolution of the field operator of noncon 
densed particles, 



With the Hamiltonian (58), we get 



(69) 



The corresponding variational derivatives give for the second-order term 



S4{v,t) 



+ y $(r - r') [^^(r')^o(r')V^i(r) + t/^^(r')V^i(r')^o(r) + ^I(r')^o(r')V'o(r)] dr' , 
for the third-order term, 

and for the fourth-order term. 



^(3) . 

J-— = / $(r - r') [4{r')M^')M^) + 4{r')M^')M^) + 4(r')V^i(r')V^i(r) dr' , 



= U{r- r')4{v')M^')M^) dr' . 
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Taking the coherent average of Eq. (69), we obtain 



2m 



+ C/-AtU'i(r)+ 



+ / $(r-r') [|,7(rO|^V'i(r) + 77*(rOV'i(rO^(r)+V'I(r077(r')77(r)+X(r,rO 
The latter equation can also be derived directly from the equation 



dt 



dr' . (70) 



(71) 



The final equation for the condensate wave function follows from Eq. (66) after averaging 
it over the space J-i. In doing this, we define for the noncondensed particles the normal 
density matrix 

pi(r,r') = <^I(r')t^i(r)>^, , (72) 

the anomalous density matrix 



(7i(r,r') = < V'i(r')V'i(r) , 
and the related diagonal densities 

pi(r) = pi(r, r) , ai{r) = ai{r, r) . 

The total particle density is 



(73) 

(74) 
(75) 



p(r) = |r7(r)|2 + pi(r) . 
Averaging Eq. (67) gives 

< X{r,r') >^,= pi(r')7?(r) + pi(r,r')7?(r') + ai(r,r')r?*(r')+ < V'I(r')V'i(r')V'i(r) ■ (76) 
Then from Eq. (66) we obtain 



' Ft ^^"^ 



2^ + ^-^'-^"^^'^+ 



+ / $(r - r') p(r')?7(r) + pi(r, r')r]{r') + ai(r, r')r]*{r')+ < '0|(r')^/'i(r')^/'i(r) dr' . 



(77) 

This is an exact equation valid for any Bose-condensed system, equilibrium or nonequilib- 
rium, uniform or nonuniform. 

For an equilibrium system, the coherent field ri{r) may be treated as not depending on 
time. Generally speaking, a stationary solution could depend on time as ?7(r, t) oc exp(— i/igt). 
This, however, would simply lead to a redefinition of the Lagrange multiplier /iq. Hence, an 
equilibrium condensate function can be defined as time independent. In that case, Eq. (77) 
yields 

" 1 

— + U{r) r]{r)+ 



sr](r) = 



2m 
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+ 



I $(r - r') [p(r')77(r) + p^{r, r')rj{r') + a^{r, r')rj*{r')+ < Vl(r')V'i(r')V'i(r) >^,] dr' , 



where the notation 



(78) 

e = /jLo + pi (79) 

is introduced, playing the role of the condensate energy per particle. Equation (78) is a 
generalized eigenproblem, which can possess the whole spectrum of the energies e and the 
related eigenfunctions ?7(r). In the case of a confining potential U{r), representing an atomic 
trap, the spectrum e is discrete, with the solutions 77 (r) representing topological coherent 
modes [38-41]. In equilibrium, the lowest-energy mode corresponds to the usual condensate. 
But in a nonequilibrium system, nonground-state condensates can be generated [38-41]. 

In a uniform system, where U{r) = 0, we have 7/(r) = rj, p(r) = p, with pi(r, r') and 
(Ti(r, r') depending on the difference r — r'. Then Eq. (78) reduces to 



erj = p^oT] + J $(r) [pi(r, 0)r] + ai{r, 0)r]*+ < V'I(0)V'i(0)V'i(r) >:f. 



dr . 



(80) 



Recall that Eqs. (77), (78), and (80) for the corresponding systems are exact, with no 
approximations involved. 



6 Dilute gas 

To illustrate the above equations for some particular cases, let us resort to the dilute-gas 
approximation, when the inequality 

Ipaf I « 1 (81) 

is valid, where p = N/V is the average density and is the s-wave scattering length. Then 
the interaction potential can be modelled by the contact form 



$(r) = $o<^(r) , 



(82) 



with 



$0 = / $(r) 
The condensate-function equation (77) becomes 



dr = An — . 
m 



2m 



+ U - e]r}{r)+ 



+$0 ip(r)r^(r) + pi(r)7^(r) + ai{r)r)*{r)+ < V'I(r)V'i(r)V'i(r) > 



(83) 



The equation of motion (70) for the field operator of noncondensed particles takes the form 
i ^ Vi(r) = (- ^ + C/ - V'i(r) + $0 [2\v{r)\'M^) + ri'{M{r) + X{r, r)] , (84) 
where 

X{r,r) = 2V^I(r)V^i(r)r;(r) + r;*(r)V^i(r)V^i(r) + V^I(r)V^i(r)V^i(r) . 
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In Eq. (83) and in what follows, we shall not specify explicitly the Fock spaces over which 
the averages are to be taken. This makes it possible to simplify the notation and cannot lead 
to confusion if, by definition, we accept the rule that each average of operators is calculated 
in that Fock space, where these operators are defined. 

Further simphfication can be done by employing the Hartree-Fock-Bogolubov (HFB) 
approximation. Then 

< V'I(r)V'i(r)Vi(r) > = . 
And the condensate-function equation (83) reduces to 



— + U-e\ 7]{r) + $0 [p{rHr) + Pi{rHr) + ai{r)r)* {r)] . (85) 



If one considers the case, when all particles are condensed, that is, p(r) 
Pi(r) = cri(r) = 0, then Eq. (85) simplifies to the equation 



l^(r 



P and 



2m 



+ C/-£ + $o|?7(r) 



77(r) . 



The latter is what is generally known as the nonlinear Schrodinger equation and in condensed- 
matter literature as the Gross-Pitaevskii equation. 

In the HFB approximation, the terms entering X(r, r) linearize as 



V'I(r)V'i(r) = V'i(r) < V'i(r) > + < i^ii^) > V'i(r)- < V'i(r) >< Vi(r) > = 



and 



4{r)Mr)Mr) = 2pi(r)V'i(r) + ai(r)V'I(r) , 
where the conservation condition (39) is taken into account. The linearized Eq. (84) is 



i|^i(r) 



2m 



+ 



U-f,] V'i(r) + $0 [2p(r)V'i(r) + <j{r)^l{r) 



where 

o-(r) = r]^{r) + (Ti(r) . 
For an equilibrium system, the condensate-function equation (85) becomes 



£7] 



- — + U(r) 
2m ^ ' 



r]{r) + $0 [p(r)?7(r) + pi(r)?7(r) + cri(r)?7*(r)] 



(86) 
(87) 

(88) 



In the case of a uniform system, when U (r) = 0, the condensate amplitude rj can be chosen 
real. From Eq. (88), taking into consideration that pi(r) = pi and (Ti(r) = Ui, we have 



£ = (P + Pl +cri)$0 • 



(89) 



Let us consider in more detail an equilibrium uniform system. The field operator of 
noncondensed particles can be expanded over the basis of plane waves 



„ikr 



16 



as 



(90) 



Because of the global symmetry of the Hamiltonian with respect to the transformation 

77(r) 77(r)e^" , V'i(r) V'i(r)e^" , (91) 
in which a is a real number, the condensate amplitude ?7(r) =7] can be taken real, 

(92) 



Po 



_ No 



The field operators and a], in the momentum representation are defined on the Fock 
space J-'{ak) generated by aj,. An important role is played by the following averages. The 
momentum distribution of particles 



Uk = < aluk > 
is called the normal average, as opposed to 



(93) 



(94) 



termed the anomalous average. The latter gives the amplitude of the process, when two 
particles are annihilated from the thermal cloud of noncondensed particles. Respectively, 
cr^ =< Ofcol^fc > is the amplitude for the process of creation of two noncondensed particles. 
With (Tfe being the amplitude of these pair processes of creation and annihilation, the absolute 
value \ak\ describes the density of binary correlated particles. 

The averages (93) and (94) define the density of noncondensed particles. 



Pi = < ^I(r)^i(r) > = 

and the amplitude of binary correlated particles, 

(7i = < Vi(r)V'i(r) > = 

The total number of noncondensed particles is 
and the number of binary correlated particles is 



1 

V 



1 



(95) 



V 



(96) 



(97) 



5i 



(98) 



So that the number of correlated pairs is Bi/2. 
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Substituting the field-operator expansion (90) into the Hamiltonian (58), wc take the 
coherent average, as a result of which the condensate operators ipo{r) are replaced by the 
coherent field r]{r). For the considered uniform system, we find the zero-order term 

= Po$o - iVo . (99) 
The first-order term H^^^ is exactly zero. The second-order term is 

The sum of Eqs. (99) and (100) corresponds to the Hamiltonian in the Bogolubov approx- 
imation [7-10]. For the higher-order terms, we obtain the following. The third-order term 
is 

H^^^ = \ly *oE' [aW+pd-p + «-p4+p«fc) , (101) 

k,p 

where the prime at the summation sign indicates that no operator in the summation corre- 
sponds to oo, so that here 

k^O, PT^O, k + p7^0. 
For the fourth-order term, we get 

H^'^ = I? E E'«W«fc+P«.-^ . (102) 

k P,q 

where the prime again means the same as earlier, that is, here it implies that 

p^O, q^O, k + p^O, q-k^O. 

The higher-order terms (101) and (102) can be treated in the HFB approximation. In so 
doing, we keep in mind that, as a consequence of the conservation condition (39), we have 

< Ofe > = (k ^ 0) . (103) 

Therefore Eq. (101) vanishes, 

H^^^ = . (104) 

And Eq. (102) reduces to 

H^^^ = $0 E hpA^k + I (TiaW-, + \ ala^kO^^ - \ {2p\ + l^il^) ^^V . (105) 
it^o ^ ^ ^ / ^ 

In this way, instead of the Hamiltonian (58), we have in the HFB approximation 

Hhfb = H^'^ + + ^(4) ^ ^106) 



18 



where is given by Eq. (99); H^^\ by Eq. (100); and H^^\ by Eq. (105). To write down 
the Hamiltonian (106) in a compact form, let us introduce the notation 

cufc = — + 2p$o - (107) 
Zm 

and, assuming that Ui is real, we define 

A = (po + ai)$o • (108) 

Then Hamiltonian (106) writes as 

Hhfb = Ehfb + ^kdW + ^ ^ (4«-fc + , (109) 

with the HFB nonoperator term 

Ehfb = H^""^ - ^ (2pi + aj) ^oV . 

Invoking the Bogolubov canonical transformation (21), discussed in Section 3, Hamilto- 
nian (109) is diagonalized to the Bogolubov form 



with the ground-state energy 



Hs^Eo + Y^kblbk, (110) 



Eq = Ehfb + XK^fe ~ '^k) 
^ k^o 



and with the spectrum 

£, = V^i-A2. (Ill) 

Separation of Bose-Einstein condensate by means of the Bogolubov shift is meaningful 
only when the momentum distribution of particles displays a singularity at the point k — 0, 
which happens if the quasiparticle spectrum, being positive, touches zero at zero momentum 
[42] . That is, the necessary condition for the occurrence of Bose-Einstein condensation is 

limefc = 0, £jk>0. (112) 

fc— >o 

This defines the chemical potential 

At = 2p$o - A , 

which, with Eq. (108), becomes 

(p + Pl -(7i)$o ■ (113) 

Combining this with Eqs. (79) and (89), we get 

IXq = £- ix = 2cTi$o • (114) 
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Substituting Eq. (113) into Eq. (107), we have 



2m 



+ (po + (7i)$o , 



which, taking into consideration Eqs 



(108) and (111), can be written as 



2m 



The quasiparticle spectrum (111) can be represented in the standard form of the Bogolubov 
spectrum 



As is evident, spectrum (115) is gapless, in agreement with the theorems by Bogolubov 
[9,10] and Hugenholtz and Pines [14]. At the same time, we have encountered no contra- 
diction with the equations of motion. Recall that the latter were derived from one given 
Hamiltonian, because of which these equations automatically preserve all conservation laws. 

The diagonal Hamiltonian (110) makes it easy to calculate different averages in the space 
J-{hk) ■ Thus, the momentum distribution of quasiparticles is 




(115) 



with the sound velocity given by the equation 




(116) 




(117) 



Taking into account that 




we find the momentum distribution of particles (93) as 




(118) 



and the anomalous average (94) as 




(119) 



For the density of noncondensed particles (95), we obtain 




(120) 
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and for the anomalous amplitude (98), we find 

In agreement with the earlier assumption, the amplitude ai can be chosen to be real. How- 
ever, the integral in Eq. (121) displays an ultraviolet divergence caused by the dilute- gas ap- 
proximation with the contact potential (82). This artificial divergence can be easily removed 
by accepting a more realistic interaction potential $(r) possessing a nontrivial momentum 
dependence of its Fourier transform One often takes the potential $(r) in the Gaussian 
form [43-45], for which $fc exponentially decreases as A; — > oo. Replacing in Eq. (108) $o 
by $fc, we get 

Ak = {po + ai)^k. (122) 
Now, the integral in Eq. (121) becomes 

„^ /llcothf^i) (123) 

The latter is convergent as soon as $fc tends to zero faster than l//c^. 

Another possibility of removing the ultraviolet divergence could be by subtracting the 
divergent term in Eq. (121), as is often done when considering hard-core particles. In such 
a manner, integral (123) could be defined as 

/•/$o "^*o\ ^.(I3^k\ dk 
a — { - — — —nr- coth 



2£fc J \ 2 J (27r)3 ■ 

Note, however, that the latter expression becomes negative, since Ek > k"^ /2m, while Eq. 
(123) is positive. 

Using notation (123), we obtain, instead of Eq. (121), the relation 

c^i = -(Po + cri)oi . 

Solving the latter, we find 

= - -f^ . (124) 

1 + a 

Because integral (123) is positive, the anomalous amplitude (124) is negative. From Eq. 
(122), we have 

A. = , (125) 

which is positive for all k. As is seen from Eqs. (124) and (125), in vicinity of the condensa- 
tion temperature T^, defined by the equation pi(Tc) = p, one has po ^ as T ^ Tc- Hence 
both (Ti and tend to zero when T —>■ Tc. But for lower temperatures T ^ Tc, where po 
can be close to p, the anomalous amplitude ai becomes of the order of pi and can even be 
much larger than the latter. Therefore the anomalous averages cannot be neglected at low 
temperatures. 

It is necessary to stress the difference between the theory, developed in this paper, and 
the standard approach to Bose-condensed systems. The principal difference is that in the 
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present theory two Lagrange multipliers, fio and /i, are introduced in Eq. (57) in order to 
preserve two normalization conditions (45) and (55). While in the standard theory, solely the 
chemical potential /i is introduced, because of this, the HFB approximation in the standard 
theory acquires an unphysical gap in the spectrum of collective excitations. If one removes 
the gap by some additional tricks, the theory becomes nonconserving. Thus, the standard 
approach always suffers from one of the defficiencies, being either gapeful or nonconserving 
[15]. 

Sometimes, to remove the gap in the collective spectrum, one sets to zero the anomalous 
average (119), ascribing this trick to Popov. First of all, such a trick is principally incorrect 
at low temperatures, when the anomalous averages become of the order or larger than the 
normal averages [46]. And, moreover, as is easy to infer from the original works by Popov 
[16-19], he has never suggested to use such an incorrect trick at low temperatures. 

In the theory, presented in this paper, there is no need to invoke any unjustified tricks. 
All averages, normal as well as anomalous, are treated on an equal footing. And the theory 
always remains both gapless and conserving. 

In conclusion to this section, it is useful to analyse the stability of matter in the HFB 
approximation. As is well known [47], the ideal uniform Bose-condensed gas is unstable. This 
instability stems from the divergent compressibility, which is directly related to anomalous 
particle fluctuations [5]. 

There have been published quite a number of papers claiming that particle fluctuations 
in the interacting Bose gas remain anomalous, being proportional to N'^^^ in the Bogolubov 
approximation. This type of anomalous behaviour does not depend on the ensemble em- 
ployed, being the same in the grand canonical and canonical Gibbs ensembles, as well as in 
the microcanonical ensemble [48]. If such anomalous fluctuations would really be present, 
this would mean the instability of the system [5,49,50]. 

However, as has been thoroughly explained in Ref. [5,49,50], there are no anomalous 
fluctuations in the Bogolubov approximation. Particle fluctuations in the Bogolubov theory 
are normal, being proportional to N. The anomalies in the fluctuations arise solely to the 
unjustifled usage of an approximation outside of the region of its validity [49,50]. 

Let us consider particle fluctuations in the HFB approximation. These fluctuations are 
characterized by the dispersion 

A2(iV) = <N'^ > - <N >^ . (126) 

A convenient expression for this dispersion is given [5,49] by the equation 

A\N)^n\i + pJ [g{r) - 1] drj , (127) 

in which the pair correlation function 

^(r) = 1 < i/;t(r)V;^(0)^(0)Vi(r) > (128) 
P 

is deflned through the fleld operators (42). 

Calculating function (128), we should keep in mind that the HFB approximation is 
equivalent to an effective mean-fleld theory, with a quadratic Hamiltonian in terms of the 
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field operators Therefore, substituting the representation (42) into Eq. (128), we 

should retain there only the terms up to the second order with respect to "01 (r). Then we 
have 

^(r) = l + ^[p,(r,0) + ai(r,0)]. 
P 

Because of the same reason, we should set here po/p — 1, since po = p — pi. Thus, we get 

Using the expressions for Uk and ak, obtained above, we find 

TN 

A^N) = -- . (130) 

That is. particle fluctuations in the HFB approximation are normal, in the same way as in 
the Bogolubov approximation [5,49,50]. 



7 Nonuniform matter 

The consideration of the previous Section 6 can be generalized to the case of nonuniform 
Bose systems with arbitrary interaction potentials $(r). As is shown in Section 5, the final 
equations of motion can be derived by projecting the Hamiltonian (58) on the space Ti. 
This implies that the coherent average of Eq. (58) is taken in the space jFg, which results 
in the replacement of the condensate field operators ipo{r) by the condensate wave function 
r){r). Keeping this projection in mind, we shall work now directly in the space J^i. 

For concreteness, we shall employ the HFB approximation, which for the fourth-order 
product of the noncondensed-particle field operators writes as 

4{r)4{r')Mr')M^) = Pi{r)4{r')Mr') + pi{r')4{r)M^) + 

+pi(r', r)V'I(r')V'i(r) + Pi(r, r')^l(v)Mr')+ 

+ai(r,r')Vl(r)V'I(r') + ar(r',r)V'i(r')Vi(r) -Pi(r)pi(r') - |pi(r,r')P - ki(r,r')P , (131) 

where the notation from Eqs. (72) to (74) is used. Then, instead of Eq. (109), we come to 
the HFB Hamiltonian for nonuniform matter 

Hhfb = Ehfb + 1 (- ^ + ^ - V'i(r) dr + J $(r - r') [pi(rOV'I(r)Vi(r)+ 

+p(r', r)Vl(r')Vi(r) + ^ <j(r, r')Vl(r')V'I(r) + ^ a*{r, r')Vi(r')V'i(r)] drdr' , (132) 
in which the notation 

p(r,r') =r?(r)7?*(r')+pi(r,r') (133) 

and 

(T(r, r') = vi^Hr') + ai(r, r') (134) 
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is introduced and where 

Ehfb = H("^- \ |$(r-r')[pi(r)pi(r') + |pi(r,r')P + |ai(r,r')r] dvdv' , (135) 



with given in Eq. (99). 

The quadratic form (132) can be diagonalized by means of the general canonical trans- 
formations [20]. For this purpose, we may expand the operators V'i(r) and "01 (i") ^ 

V'i(r) = EhW^'fc + ^:(r)fcI] , 

^\{v) = Y.{ul{r)h\ + v^{v)h^ , (136) 

k 

where k is a set of quantum numbers. The transformations, inverse to Eq. (136), are 



«fe(r)^i(r) - vl{r)ijl{r)] dr , 



bl 



Mfc(r)^I(r) - Vfc(r)^i(r) 



dr . 



(137) 



The operators V'i(r) and bk are assumed to satisfy the Bose commutation relations, which 
imposes on the coefficient functions Uk{r) and Vk{r) the following restrictions: 

EiMrKir') - vl{r)u,{r')] = , 

k 

Y,[u,{r)ul{r')-vl{r)v,{r')] = 6{r-r'), 



[uk{r)vp{r) - Vk{r)up{r)] dr = , 
j K{r)up{r) - vl{r)vp{r)\ dr = 5kp ■ 



(138) 



Let us define 
uir, r') ^ 

and 



2m 



+ U{r) - /X + / $(r - r')p(r') dr' 



5{r - r') + $(r - r')p(r, r') (139) 

(140) 



A(r,r') = $(r-r')a(r,r') . 
Then the diagonalization condition can be represented as a system of equations 



j [uj{r, r')uk{r') + A(r, r>fc(r')] dr' = ekUk{r) , 
j [a;*(r, r')vk{r') + A*(r, r')uk{r')\ dr' = -ekVk{v) , 



(141) 



which is a variant of the Bogolubov - de Gennes equations. The HFB Hamiltonian (132) 
acquires the diagonal form 

^B = £^o + E^A, (142) 

k 
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which is analogous to Eq. (110), with the difference that here k is not momentum but a set 
of quantum numbers labelling the functions Uk{r) and Vk{r). And the ground-state energy 
is 

Eo = Ehfb -Y.^kj \vk{^)? dv . (143) 

k 

Calculating the averages in the space T{bk), with the Hamiltonian (142), we find the 
normal density matrix 

pi(r, r') = Yl [^kUk{r)ul{r') + (1 + nk)vl{r)vk{r')] (144) 

k 

and the anomalous density matrix 

(7i(r, r') = ^ [KkUk{r)vl{r') + (1 + nk)vl{r)uk{r')] , (145) 

k 

defined in Eqs. (72) and (73), with tt^ from Eq. (117). The corresponding diagonal elements 
(74) give the density of noncondensed particles 

Pi(r) = E hl«fc(r)|' + (1 + ^k)\vk{r)f] (146) 

k 

and the anomalous amplitude 

a,{r)=J2{l + 27Tk)uk{r)vl{r) . (147) 

k 

These general equations are valid for an arbitrary external potential U (r) and for any pair 
interaction potential $(r). 

One can always return to the uniform system by taking 

Uk(r) = Ukipk(r) , Vk(r) = Vkipk{r) , (148) 

where (pk{'^) is a plane wave. By expanding in the plane waves the kernel (139) as 

a;(r,r')=E^fc^fc(r)^:(r') (149) 

k 

and Eq. (140) as 

A{r,r') = J2AkMrM{r'): (150) 

k 

we obtain cuk and A^, which are equivalent to Eqs. (107) and (108), respectively. 

In this way, any system with Bose-Einstein condensate can be described by a self- 
consistent theory, which is both conserving and gapless. The basis for this theory is an 
accurate usage of the nonequivalent operator representations, associated with canonical com- 
mutation relations, and the derivation of equations of motion taking into consideration all 
conditions characterizing the Bose-condensed system. In particular, taking into account two 
normalization conditions (45) and (55) requires to introduce two Lagrange multipliers //q 
and /X. 
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